The attitude of a space vehicle is defined by the relationship between its body reference frame and an inertial reference frame. There are several parameter sets that can be used to define the relationship between the two coordinate frames. The most popular are the Euler rotation vector, the direction cosine transformation matrix (DCM), Euler angles, and quaternions. The rate of attitude change is defined by the angular rate vector. Each attitude parameter set can be propagated when the angular rate is known. Thus, there is an associated differential equation driven by the angular rate vector for each parameter set.
Attitude propagation is accurate and easy, if the angular rate vector points in a fixed direction in space. In this special case, the angular rate vector can be integrated to obtain the Euler rotation vector. The other attitude parameters can be easily obtained from the Euler rotation vector. However, most attitude propagation cases of interest involve angular rate vectors that change direction as a function of time. In these cases, the angular rate vector cannot be integrated to yield the Euler rotation vector. One could assume that the angular rate vector remains constant over a short time interval and define a series of incremental angles. The incremental rotation angles could then be used to update the attitude. The smaller the time increment, the smaller the error over a fixed time duration. However, there will always be an error that grows with time. This error is often referred to as “coning error” because it appears most clearly when the angular rate vector executes a coning motion.
Propagating the attitude using the Euler rotation vector is popular because it contains only three parameters and does not have to be renormalized periodically as do quaternions or DCMs. In addition, the differential equation for the Euler rotation vector has a coning term. It is known in the prior art to evaluate numerically the coning term to reduce or eliminate coning error, which is important for vehicles having strapdown Inertial Navigation Units (INUs). These vehicles require that the attitude propagation algorithms be fast and accurate, to provide the correct attitude in real time or faster. Once the Euler rotation vector is propagated, it can be used to create the associated DCM or quaternion as needed.
As examples of prior art techniques for propagating attitude, the attitude of a vehicle can be specified by a transformation matrix {right arrow over (U)} from the vehicle's coordinate frame to an inertial reference coordinate. If {right arrow over (X)}B is a vector in the vehicle body, it transforms to the inertial frame by the transformation,{right arrow over (X)}I={right arrow over (U)}{right arrow over (X)}B  (1)where {right arrow over (U)} is the direction cosine transformation matrix and {right arrow over (X)}I is the vector in the inertial frame. The direction cosine transformation matrix {right arrow over (U)} can be expressed as:
                              U          ⟶                =                  (                                                                      1                  +                                                            (                                                                        e                          3                          2                                                +                                                  e                          2                          2                                                                    )                                        ⁢                    C                                                                                                -                                      (                                                                                            e                          1                                                ⁢                                                  e                          2                                                ⁢                        C                                            +                                                                        e                          3                                                ⁢                        S                                                              )                                                                                                                                          e                      2                                        ⁢                    S                                    -                                                            e                      1                                        ⁢                                          e                      3                                        ⁢                    C                                                                                                                                                                  e                      3                                        ⁢                    S                                    -                                                            e                      1                                        ⁢                                          e                      2                                        ⁢                    C                                                                                                1                  +                                                            (                                                                        e                          1                          2                                                +                                                  e                          3                          2                                                                    )                                        ⁢                    C                                                                                                -                                      (                                                                                            e                          3                                                ⁢                                                  e                          2                                                ⁢                        C                                            +                                                                        e                          1                                                ⁢                        S                                                              )                                                                                                                        -                                      (                                                                                            e                          3                                                ⁢                                                  e                          1                                                ⁢                        C                                            +                                                                        e                          2                                                ⁢                        S                                                              )                                                                                                                                          e                      1                                        ⁢                    S                                    -                                                            e                      2                                        ⁢                                          e                      3                                        ⁢                    C                                                                                                1                  +                                                            (                                                                        e                          1                          2                                                +                                                  e                          2                          2                                                                    )                                        ⁢                    C                                                                                )                                    (        2        )            where {right arrow over (e)} is the unit vector in the direction of an axis of rotation known as the Euler axis, θ is the angle of rotation about the Euler axis, C=cos(θ)−1, and S=sin(θ).
Euler's' rotation theorem, which is graphically illustrated in FIG. 1, states that any transformation between two coordinate frames can be represented by a single rotation of angle θ about the Euler rotation axis. Thus, an Euler rotation vector, {right arrow over (θ)}, can be defined to represent an arbitrary transformation where:{right arrow over (θ)}=θx{right arrow over (e)}1+θy{right arrow over (e)}2+θz{right arrow over (e)}3  (3)Given the Euler rotation vector, {right arrow over (U)} can be obtained using equation (2) above.
The direction cosine transformation matrix can also be expressed in terms of the quaternion elements as follows:
                              U          ⟶                =                  (                                                                                          q                    1                    2                                    -                                      q                    2                    2                                    -                                      q                    3                    2                                    +                                      q                    4                    2                                                                                                2                  ⁢                                      (                                                                                            q                          1                                                ⁢                                                  q                          2                                                                    -                                                                        q                          3                                                ⁢                                                  q                          4                                                                                      )                                                                                                2                  ⁢                                      (                                                                                            q                          1                                                ⁢                                                  q                          3                                                                    +                                                                        q                          2                                                ⁢                                                  q                          4                                                                                      )                                                                                                                        2                  ⁢                                      (                                                                                            q                          1                                                ⁢                                                  q                          2                                                                    +                                                                        q                          3                                                ⁢                                                  q                          4                                                                                      )                                                                                                                    -                                          q                      1                      2                                                        +                                      q                    2                    2                                    -                                      q                    3                    2                                    +                                      q                    4                    2                                                                                                2                  ⁢                                      (                                                                                            q                          2                                                ⁢                                                  q                          3                                                                    -                                                                        q                          1                                                ⁢                                                  q                          4                                                                                      )                                                                                                                        2                  ⁢                                      (                                                                                            q                          1                                                ⁢                                                  q                          3                                                                    -                                                                        q                          2                                                ⁢                                                  q                          4                                                                                      )                                                                                                2                  ⁢                                      (                                                                                            q                          2                                                ⁢                                                  q                          3                                                                    +                                                                        q                          1                                                ⁢                                                  q                          4                                                                                      )                                                                                                                    -                                          q                      1                      2                                                        -                                      q                    2                    2                                    +                                      q                    3                    2                                    +                                      q                    4                    2                                                                                )                                    (        4        )            The first three quaternion elements form the vector component of the quaternion. This vector component is aligned with the Euler rotation vector and the fourth quaternion component is a scalar. The quaternion elements are defined by equations (5) to (9) below:
                                          Q            ⟶                    ⁡                      (            θ            )                          =                              (                                                                                q                    1                                                                                                                    q                    2                                                                                                                    q                    3                                                                                                                    q                    4                                                                        )                    =                      (                                                                                q                    ⟶                                                                                                                    q                    4                                                                        )                                              (        5        )                                          q          1                =                              e            1                    ⁢                      sin            ⁡                          (                              θ                2                            )                                                          (        6        )                                          q          2                =                              e            2                    ⁢                      sin            ⁡                          (                              θ                2                            )                                                          (        7        )                                          q          3                =                              e            3                    ⁢                      sin            ⁡                          (                              θ                2                            )                                                          (        8        )                                          q          4                =                  cos          ⁡                      (                          θ              2                        )                                              (        9        )            
An attitude transformation can also be represented by a sequence of rotations about the axes of a coordinate system. A typical Euler rotation sequence is yaw, pitch, roll, where yaw is about the z-axis, pitch is about the y-axis, and roll is about the x-axis. The total transformation is the product of the individual transformation matrices and is given by:{right arrow over (U)}={right arrow over (U)}({right arrow over (θ)}Z){right arrow over (U)}({right arrow over (θ)}y){right arrow over (U)}({right arrow over (θ)}x)  (10)The Euler rotation sequences, the Euler rotation vectors, and the quaternions can all be transformed into their corresponding DCMs.
The vehicle attitude can be propagated when the angular rates are known. The first order time dependent differential equations for quaternion parameters are given by:
                              (                                    ⅆ                              Q                ⟶                                                    ⅆ              t                                )                =                              (                          1              2                        )                    ⁢                      (                                                            0                                                  r                                                                      -                    q                                                                    p                                                                                                  -                    r                                                                    0                                                  p                                                  q                                                                              q                                                                      -                    p                                                                    0                                                  r                                                                                                  -                    p                                                                                        -                    q                                                                                        -                    r                                                                    0                                                      )                    ⁢                      Q            ⟶                                              (        11        )            where p, q, r are roll, pitch, and yaw rates about the x, y, and z vehicle axes, respectively. The elements of the direction cosine matrix can be updated based on the set of differential equations given by:
                              (                                    ⅆ                              U                ⇀                                                    ⅆ              t                                )                =                              U            ⟶                    ⁡                      (                                                            0                                                                      -                    r                                                                    q                                                                              r                                                  0                                                                      -                    p                                                                                                                    -                    q                                                                    p                                                  0                                                      )                                              (        12        )            
The quaternions can be integrated using equation (11) and used in equation (4) to obtain the direction cosine transformation matrix, which can be used to transform vectors between vehicle body and inertial frames. Alternatively, the elements of the direction cosine matrix can be integrated directly using equation (12).
Another method to propagate vehicle attitude is by combining a series of small rotations associated with a set of small time steps. It can be assumed that the angular rate vector remains fixed over the small time intervals. Each rotation results in an incremental rotation about the angular rate axis, which is considered an Euler rotation axis. The incremental rotation vector and transformation matrix are given by:
                              Δ          ⁢                      θ            ⟶                          =                              (                                                                                Δθ                    x                                                                                                                    Δθ                    y                                                                                                                    Δθ                    z                                                                        )                    =                      (                                                                                p                    ⁢                                                                                  ⁢                    Δ                    ⁢                                                                                  ⁢                    t                                                                                                                    q                    ⁢                                                                                  ⁢                    Δ                    ⁢                                                                                  ⁢                    t                                                                                                                    r                    ⁢                                                                                  ⁢                    Δ                    ⁢                                                                                  ⁢                    t                                                                        )                                              (        13        )                                V        =                  (                                                    0                                                                                  -                    Δ                                    ⁢                                                                          ⁢                                      θ                    z                                                                                                Δθ                  y                                                                                                      Δθ                  z                                                            0                                                              -                                      Δθ                    x                                                                                                                        -                                      Δθ                    y                                                                                                Δθ                  x                                                            0                                              )                                    (        14        )            The updated transformation matrix is given by:{right arrow over (U)}(t+Δt)={right arrow over (U)}(t)+{right arrow over (U)}(t){right arrow over (V)}  (15)The associated quaternion is given by:
                                          Q            ⇀                    ⁡                      (                          t              +                              Δ                ⁢                                                                  ⁢                t                                      )                          =                                            Q              ⟶                        ⁡                          (              t              )                                +                                    (                              1                2                            )                        ⁢                          (                                                                    0                                                                              Δθ                      z                                                                                                                          -                        Δ                                            ⁢                                                                                          ⁢                                              θ                        y                                                                                                                        Δθ                      x                                                                                                                                  -                                              Δθ                        z                                                                                                  0                                                                              Δθ                      x                                                                                                  Δθ                      y                                                                                                                                  Δθ                      y                                                                                                  -                                              Δθ                        x                                                                                                  0                                                                              Δθ                      z                                                                                                                                  -                                              Δθ                        x                                                                                                                        -                                              Δθ                        y                                                                                                                        -                                              Δθ                        z                                                                                                  0                                                              )                        ⁢                                          Q                ⇀                            ⁡                              (                t                )                                                                        (        16        )            
This attitude propagation procedure has been enhanced to improve attitude propagation accuracy. With reference to FIG. 5, the first order differential equation for the Euler rotation vector, {right arrow over (θ)}, involves the angular rate vector, {right arrow over (ω)}, and is given by:
                              θ                      ⟶            .                          =                              ω            ⟶                    +                                    1              2                        ⁢                          θ              ⟶                        ×                          ω              ⟶                                +                                    (                              1                -                                                      θ                    ⁢                                                                                  ⁢                                          sin                      ⁡                                              (                        θ                        )                                                                                                  2                    ⁢                                          (                                              1                        -                                                  cos                          ⁡                                                      (                            θ                            )                                                                                                                                                          )                        ⁢                          θ              ⟶                        ×                          (                                                θ                  ⟶                                ×                                  ω                  ⟶                                            )                                                          (        17        )            
Equation (13) is equivalent to equation (17) integrated over a short time interval without the last two terms on the right hand side in equation (17). The last two terms on the right hand side in equation (17) account for the change direction of the angular rate vector. In many implementations, the last term on the right hand side in equation (17) is omitted. The prior art comprises numerical methods to integrate equation (17). These methods are often referred to as “coning correction algorithms.” Despite the improved accuracy that is obtained when using coning correction algorithms, there are still errors associated with coning correction algorithms.